<html>
  <head>
    <meta http-equiv="Content-Type" content="text/html; charset=UTF-8">
    <title>phc</title>
  </head>
  <body bgcolor="#FFFFFF">
    <center>Scilab Function</center>
    <div align="right">Last update : April 1993</div>
    <p>
      <b>phc</b> -  Markovian representation</p>
    <h3>
      <font color="blue">Calling Sequence</font>
    </h3>
    <dl>
      <dd>
        <tt>[H,F,G]=phc(hk,d,r)  </tt>
      </dd>
    </dl>
    <h3>
      <font color="blue">Parameters</font>
    </h3>
    <ul>
      <li>
        <tt>
          <b>hk</b>
        </tt>: hankel matrix</li>
      <li>
        <tt>
          <b>d</b>
        </tt>: dimension of the observation</li>
      <li>
        <tt>
          <b>r</b>
        </tt>: desired dimension of the state vector for the approximated model</li>
      <li>
        <tt>
          <b>H, F, G</b>
        </tt>: relevant matrices of the Markovian model</li>
    </ul>
    <h3>
      <font color="blue">Description</font>
    </h3>
    <p>
    Function which computes the matrices <tt>
        <b>H, F, G</b>
      </tt> of a Markovian 
    representation by the principal hankel
    component approximation method, from the hankel matrix built
    from the covariance sequence of a stochastic process.</p>
    <h3>
      <font color="blue">Examples</font>
    </h3>
    <pre>

//
//This example may usefully be compared with the results from 
//the 'levin' macro (see the corresponding help and example)
//
//We consider the process defined by two sinusoids (1Hz and 2 Hz) 
//in additive Gaussian noise (this is the observation); 
//the simulated process is sampled at 10 Hz.
//
t=0:.1:100;rand('normal');
y=sin(2*%pi*t)+sin(2*%pi*2*t);y=y+rand(y);plot(t,y)
//
//covariance of y
//
nlag=128;
c=corr(y,nlag);
//
//hankel matrix from the covariance sequence
//(we can choose to take more information from covariance
//by taking greater n and m; try it to compare the results !
//
n=20;m=20;
h=hank(n,m,c);
//
//compute the Markov representation (mh,mf,mg)
//We just take here a state dimension equal to 4 :
//this is the rather difficult problem of estimating the order !
//Try varying ns ! 
//(the observation dimension is here equal to one)
ns=4;
[mh,mf,mg]=phc(h,1,ns);
//
//verify that the spectrum of mf contains the 
//frequency spectrum of the observed process y
//(remember that y is sampled -in our example 
//at 10Hz (T=0.1s) so that we need 
//to retrieve the original frequencies through the log 
//and correct scaling by the frequency sampling)
//
s=spec(mf);s=log(s);
s=s/2/%pi/.1;
//
//now we get the estimated spectrum
imag(s),
//
 
  </pre>
    <h3>
      <font color="blue">See Also</font>
    </h3>
    <p>
      <a href="levin.htm">
        <tt>
          <b>levin</b>
        </tt>
      </a>,&nbsp;&nbsp;</p>
  </body>
</html>
